%matplotlib widget
import matplotlib.pyplot as plt
import numpy as np
import sympy
from src.v2.impl.conditions import StepCountCondition, PrecisionCondition, AbsolutePrecisionCondition
from src.v2.impl.methods import CoordinateDescent, GoldenRatioMethod, NewtonWolfe
from src.v2.impl.metrics import StepCount, CallCount, GradientCallCount, HessianCallCount, PrecisionCount, \
AbsolutePrecisionCount, AbsolutePrecision, MinAbsolutePrecision
from src.v2.impl.oraculs import LambdaOracul, SymbolOracul
from src.v2.runner.debug import FULL_DEBUG
from src.v2.runner.runner import Runner, FULL_VISUALIZE, NO_VISUALIZE, FULL_ANIMATION
from src.v2.visualization.animation import Animator
from src.v2.runner.runner import TABLE
from src.v2.impl.methods import GradientDescent, ScipyMethod, Newton, NewtonBase
from IPython.display import display, HTML
display(HTML("<style>pre { white-space: pre !important; }</style>"))
def print_points(data):
for i in data:
print(i[0] + i[1])
Задачей данной работы нам видится исследование эффективности методов на различных функциях с целью выявления сильных и слабых их сторон. Особое внимание, в силу тематики работы, уделяется методам Ньютона и квазиньютоновским методам
Для начала, зададим тестовые функции. Согласно рекомендации к выполнению, конечно, среди тестовых функций, будет функция Розенброка. Кроме того, для проверки корректности поведения и базовой работоспособности, добавим в подборку базовую функцию - квадратичную и не представляющую трудностей для методов Ньютона (ожидаемое поведение - схождение за 2 шага с достижением минимума на первом) и чуть более сложный случай - функцию Била. В качестве неполиномиальных функций для рассмотрения возьмём, для разминки, сумму экспонент, более же содержательным экзепляром будет функция hard_non_poly. Также стоит обратить внимание на то, что многие функции дублированы. Это связано с реализациями различных способов вычисления. Функции, помеченнные _absolute производят вычисление производных символьно посредством библиотеки sympy. Это требует несколько большего времени при вычислении из-за необходимости подстановки значений в символьное выражение, однако даёт максимально возможную точность при вычислениях. Оракулы же без данной пометки производят вычисления производной и гессиана посредством средств numdifftools разностным методом, потому дают приближё.нное значение
Точки старта будут задаваться в каждом окне, запускающем вычисление в переменной point. В анализе результатов, дабы не загромождать их лишним текстом, повторяться выбор точки не будет
Красной точкой на графиках отображена точка остановки алгоритма, анимации отключены для облегчения вопсриятия и повышения производительности. В случаях, когда визуализация включена, каждой функции соответствует 2 графика - 3d и линии уровня, расположенных в соответствии с порядком функций в таблице
x, y, z = sympy.symbols('x, y, z')
base = LambdaOracul(lambda x, y: (x - 10) ** 2 + (y - 5) ** 2)
base_absolute = SymbolOracul((x - 10) ** 2 + (y - 5) ** 2, ['x', 'y'])
bute = LambdaOracul(lambda x, y: (x + 2 * y - 7)**2 + (2 * x + y - 5)**2)
bute_absolute = SymbolOracul((x + 2 * y - 7)**2 + (2 * x + y - 5)**2, ['x', 'y'])
rosenbrok = LambdaOracul(lambda x, y: (1 - np.float64(x))**2 + 100 * (np.float64(y) - np.float64(x)**2)**2)
rosenbrok_absolute = SymbolOracul((1 - x)**2 + 100 * (y - x**2)**2, ['x', 'y'])
bil = LambdaOracul(lambda x, y: (1.5 - np.float64(x) + np.float64(x) * np.float64(y))**2 + (2.25 - np.float64(x) + np.float64(x) * np.float64(y)**2)**2 + (2.625 - np.float64(x) + np.float64(x) * np.float64(y)**3)**2)
hard_non_poly = LambdaOracul(lambda x, y: np.sqrt(np.float64(x)**2+np.float64(y)**2) - np.sin(np.float64(x) + np.float64(y))**2)
macCormic = LambdaOracul(lambda x,y: np.sin(x + y) + (np.float64(x) - np.float64(y))**2 - 1.5 * np.float64(x) + 2.5 * np.float64(y) + 1)
base_non_poly = LambdaOracul(lambda x, y: np.exp(np.abs(x, dtype=np.float64), dtype=np.float64) + np.exp(np.abs(y, dtype=np.float64), dtype=np.float64))
min_bil = [3, 0.5]
min1 = np.array([10, 5])
oraculs_rosenbrok = [rosenbrok]
min2 = [1, 1]
oraculs_bute = [bute_absolute]
min3 = [1, 3]
oraculs_hard = [base_non_poly, hard_non_poly]
min4 = [0, 0]
В качестве метрик для сравнения во всех последующих испытаниях будут использованы число шагов, совершённых алгоритмом до остановки(StepCount), число вызовов функции (CallCount), число вызовов градиента (GradientCallCount), число вызовов Гессиана (HessianCallCount), а также число шагов до достижения реальной требуемой точности (AbsolutePrecisionCount) и реальная точность на момент схождения алгоритма (AbsolutePrecision) - последние две метрики для рассчёта используют знание о реальном минимуме функции, потому стоит заметить, что метрики являются системами-наблюдателями и на работу методов не влияют.
Ещё одним важным замечанием является то, что ввиду некоторой инородности для нашей системы методов scipy, метрика AbsolutePrecisionCount работает для них не вполне корректно, выдавая в случае схождения 1. Данная особенность связана с деталями реализации. В отчёте же данная метрика оставлена в большей степени в качества маркера. Значение Undefined означает, что методу сойтись к истинному минимому не удалось
Говоря про точность, далее будут использованы 3 основном: minPrec=1e-5, maxPrec=1e-11 и defPrec=1e-7. Именно последняя будет использована в большинстве случаев. Также здесь стоит упомянуть, что название maxPrec обусловлено не столько реальными возможностями алгоритмов, сколько тем фактом, что нецелочисленные типы хранят лишь 14 значащих цифр, потому достижение большей точности может вызывать чисто технические сложности, особенно если минимумы находятся в точках, значительно больших нуля.
metrics_base = [StepCount(), CallCount(), GradientCallCount(), HessianCallCount()]
minPrec=1e-5
defPrec=1e-7
maxPrec=1e-11
abs_metrics1 = metrics_base + [AbsolutePrecisionCount(defPrec, min1), AbsolutePrecision(min1)]
abs_metrics2 = metrics_base + [AbsolutePrecisionCount(defPrec, min2), AbsolutePrecision(min2)]
abs_metrics3 = metrics_base + [AbsolutePrecisionCount(defPrec, min3), AbsolutePrecision(min3)]
abs_metrics4 = metrics_base + [AbsolutePrecisionCount(minPrec, min4), AbsolutePrecision(min4)]
abs_metr_bil = metrics_base + [AbsolutePrecisionCount(minPrec, min_bil), AbsolutePrecision(min_bil)]
animations = [Animator()]
Остановка методов происходит по достижении искомой точности (PrecisionCondition) или лимита шагов (равного, как правило, 500)
conditions = [StepCountCondition(500), PrecisionCondition(defPrec)]
modules1 = animations + abs_metrics1 + conditions
modules2 = animations + abs_metrics2 + conditions
modules3 = animations + abs_metrics3 + conditions
modules4 = animations + abs_metrics4 + [StepCountCondition(1000), PrecisionCondition(minPrec)]
modules_bil = animations + abs_metr_bil + [StepCountCondition(500), PrecisionCondition(minPrec)]
Первая серия экспериментов будет представлять из себя сравнение метода Newton-CG, представленного в библиотеке scipi, и методов Ньютона, реализованных нами (с фиксированным и плавающим шагами. В методе с плавающим шагом использован метод золотого сечения, сравнение эффективности с методом Вольфе, предлагаемым к реализации в первом дополнительным задании представлено в последующих сериях, относящихся к дополнительному заданию 1
Краткое описание методов:¶
- NewtonBase - классический метод Ньютона с фиксированным шагом. Параметризуется параметром learning_rate, отвечающим за длину шага по направлению предполагаемого минимума. Каждая следующая точка вычислется по правилу x = x_prev - hess^(-1)(x_prev)grad(x_prev)learning_rate, где grad - значение градиента функции в переданной точке, а hess - значение гессиана. В случае, если не существует hess^-1, значения на главной диагонали увеличиваются на eps=1e-7, что, как правило, приводит к выбору в качестве направления направление градиента
- Newton - метод Ньютона с перменным шагом. Параметризуется параметрами learning_rate - максимальная длина шага, method - метод, используемый для поиска минимума на луче, aprox_dec - точность приближения для внуреннегоо метода. Луч выбирается аналогичным NewtonBase образом. В тестах в качестве метода используется метод золотого сечения. Это связано с тем, что метод Newton-CG использует в свой реализации правило Армихо, правило Вольфе же вынесено в дополнительное задание 1, потому, для более полного демонстрация различий (и без того не слишком явных, как будет заметно далее), зависящих от метода одномерного поиска была выбрана данная реализация
- Newton-CG - метод Ньютона, реализованный в библиотеке scipy. Внутренняя реализация для поиска по направлению использует правило Армихо, вычисление значения, обратного Гессиану, выполняется приближённо. В качестве параметризации имеются параметры c1 - для правила Армихо и c2 - для правила Кривизны
- BFGS - квазиньютоновский метод, представленный в библиотекке scipy. Параметризация и особенность реализации аналогичны методу Newton-CG. Отличие заключается в том, что вместо вычисления Гессиана вычислется его приближённое значение с использованием градиента
newton_methods = [
NewtonBase(1),
Newton(aprox_dec=maxPrec),
ScipyMethod("Newton-CG")
]
point = np.array([1021.2, 3200.5])
result = Runner.run(newton_methods, [base_absolute], point, modules1, precision=defPrec, **TABLE, **FULL_VISUALIZE)
SymbolOracul +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Method name | StepCount | CallCount | GradientCallCount | HessianCallCount | AbsolutePrecisionCount(1e-07) | AbsolutePrecision | +=======================================+=============+=============+=====================+====================+=================================+=====================+ | NewtonBase(1) | 2 | 0 | 2 | 2 | 1 | 0 | +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Newton(3,GoldenRatioMethod,eps=1e-11) | 2 | 114 | 2 | 2 | 1 | 8.88178e-16 | +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | ScipyMethod(newton-cg) | 2 | 2 | 2 | 2 | 1 | 0 | +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+
Вычисления производились на точных значениях производной. Поведение метода совпадает с ожидаемым: все реализации сошлись за 2 шага, достигнув точного минимума на первом (В случае с реализацией метода Ньютона на основе метода золотого сечения точный минимум достигнут не был, ведь вместо перемещения по шагу происходит поиск минимума на луче -> приближённое значение. Хотя теоритически и возможно случачйное попадание в точное знаение минимума) Стоит заметить, что хотя
point = np.array([102.1, 320.5])
result = Runner.run(newton_methods, oraculs_bute, point, modules3, precision=defPrec, **TABLE, **NO_VISUALIZE)
SymbolOracul +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Method name | StepCount | CallCount | GradientCallCount | HessianCallCount | AbsolutePrecisionCount(1e-07) | AbsolutePrecision | +=======================================+=============+=============+=====================+====================+=================================+=====================+ | NewtonBase(1) | 2 | 0 | 2 | 2 | 1 | 0 | +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Newton(3,GoldenRatioMethod,eps=1e-11) | 2 | 114 | 2 | 2 | 1 | 1.11022e-16 | +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | ScipyMethod(newton-cg) | 5 | 9 | 9 | 9 | 1 | 2.22045e-16 | +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+
Функция Бута так же не вызвала затруднений у методов. Хотя стоит заметить, что метод Ньютона с золотым сечением внутри сошёлся быстрее реализации scipy
point = np.array([102.1, 320.5])
result = Runner.run(newton_methods, oraculs_rosenbrok, point, modules2, precision=defPrec, **TABLE, **FULL_VISUALIZE)
LambdaOracul +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Method name | StepCount | CallCount | GradientCallCount | HessianCallCount | AbsolutePrecisionCount(1e-07) | AbsolutePrecision | +=======================================+=============+=============+=====================+====================+=================================+=====================+ | NewtonBase(1) | 7 | 0 | 7 | 7 | 6 | 2.49247e-14 | +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Newton(3,GoldenRatioMethod,eps=1e-11) | 124 | 7068 | 124 | 124 | 123 | 7.05057e-15 | +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | ScipyMethod(newton-cg) | 75 | 276 | 276 | 149 | 1 | 7.07328e-11 | +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+
Данное тестирование проходило на функции Розенброка, производные вычислялись численными методами. Несложно заметить, что все методы сошлись до искомой точности. Как и в прошлом эксперименте, метод, делающий полный шаг сошёлся до точности более 14 знаков, притом потребовалось для этого лишь 7 шагов. Это обусловлено тем, что хотя функция и не в точности квадратичная, вид её крайне на неё похож, кроме того, функция полиномиальна. Вместе это позволяет получать подобную точность
Другие реализации метода так же ожидаемым образом достигли искомой точности. Стоит заметить, что метод Ньютона как в данном, так и в предыдущем испытании сошёлся до несколько большей точности, чем было запрошено - это обусловлено спецификой раобты метода золотого сечения и является стабильной особенностью.
Также стоит заметить, что на простых для метода функциях испольование подобного метода одномерного поиска излишне. Меньшая точность шага вынуждает делать большее число шагов, кроме того на каждом шаге запрашивается значительное число значений функции. В совокупности это ставит под сомнение вычислительную эффективность подобной реализации.
Реализация scipi, в силу своих особенностей, запрашивает меньшее число значений, однако требует несколько большего числа вычислений градиента и гессиана. Но и потери в точности могут быть весьма весомы. Ниже приведён пример, когда реализованные с честным вычислением Гессиана метода сошлись к минимуму, в то время как Newton-CG ушёл от него
point = np.array([-100, 20.5])
methods = [
NewtonBase(1),
Newton(learning_rate = 10, aprox_dec=maxPrec),
ScipyMethod("Newton-CG")
]
result = Runner.run(newton_methods, oraculs_rosenbrok, point, modules2, precision=defPrec, **TABLE, **FULL_VISUALIZE)
LambdaOracul +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Method name | StepCount | CallCount | GradientCallCount | HessianCallCount | AbsolutePrecisionCount(1e-07) | AbsolutePrecision | +=======================================+=============+=============+=====================+====================+=================================+=====================+ | NewtonBase(1) | 7 | 0 | 7 | 7 | 6 | 7.46747e-14 | +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Newton(3,GoldenRatioMethod,eps=1e-11) | 124 | 7068 | 124 | 124 | 123 | 4.36954e-15 | +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | ScipyMethod(newton-cg) | 201 | 401 | 401 | 400 | Undefined | 21.0484 | +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+
Тесты на неполиномиальных функциях. Первая из функций - сумма экспонент, вторая же будет подробнее описана при разборе результатов
point = np.array([-10, 20.5])
newton_methods = [
NewtonBase(1),
Newton(learning_rate = 5, aprox_dec=maxPrec),
ScipyMethod("Newton-CG")
]
result = Runner.run(newton_methods, oraculs_hard, point, modules4, precision=minPrec, **TABLE, **FULL_VISUALIZE)
LambdaOracul +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Method name | StepCount | CallCount | GradientCallCount | HessianCallCount | AbsolutePrecisionCount(1e-05) | AbsolutePrecision | +=======================================+=============+=============+=====================+====================+=================================+=====================+ | NewtonBase(1) | 999 | 0 | 999 | 999 | Undefined | 0.5 | +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Newton(5,GoldenRatioMethod,eps=1e-11) | 6 | 348 | 6 | 6 | 5 | 2.39462e-15 | +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | ScipyMethod(newton-cg) | 14 | 30 | 30 | 26 | 1 | 3.55239e-06 | +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ LambdaOracul +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Method name | StepCount | CallCount | GradientCallCount | HessianCallCount | AbsolutePrecisionCount(1e-05) | AbsolutePrecision | +=======================================+=============+=============+=====================+====================+=================================+=====================+ | NewtonBase(1) | 186 | 0 | 186 | 186 | Undefined | 2.2518e+15 | +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Newton(5,GoldenRatioMethod,eps=1e-11) | 4 | 232 | 4 | 4 | Undefined | 7.49736 | +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | ScipyMethod(newton-cg) | 3 | 8 | 8 | 5 | Undefined | 7.49736 | +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+
print_points(result[0])
('no_show_value', array([-3.6722549e-17, -5.0000000e-01]), 'Method name', 'NewtonBase(1)')
('no_show_value', array([-7.77100915e-16, -2.26501701e-15]), 'Method name', 'Newton(5,GoldenRatioMethod,eps=1e-11)')
('no_show_value', array([ 2.51191848e-06, -2.51191848e-06]), 'Method name', 'ScipyMethod(newton-cg)')
На неполиномиальных функциях методы расккрываются уже с другой стороны. Безоговорочный лидер предыдущих испытний - метод Ньютона с шагом 1 на функции суммы экспонент не сумел даже сойтись за отведённое число шагов. Подобное поведение во многом ожидаемо. Функция в окресности нуля крайне плоская, кроме того, в виду своей неполиномиальности, её минимум не вычисляется столь грубым приближением. Наиболее же эффективным на простой неполиномиальной функции оказался метод Ньютона, основанный на золотом сечении, сойдясь всего за 4 итерации. Стоит заметить, что 4 итерации - это всего 4 вычисления Гессиана против 29 реализации scipy. Точность же практически в 2 раза превысила запрошенную. Функция hard_non_poly оказалась показательной в том отношении, что, хотя, как показала первая функция, методы способны работать с неполиномиальными функциями (при условии их унимодальности), даже небольшая неполиномиальность, особенно порождающая локальные минимумы, может оказаться непреодолимой для методов. Как несложно заметить, данная функция представляет собой воронку с единственным глобальным минимумом в точке (1, 1). За счёт члена sin(x + y) порождаются углубления. NewtonBase сойтись не сумел вовсе, оставшиеся методы сошлись к локальному минимуму (что примечательно, одному и тому же)
Таким образом, на функциях простого вида наиболее эффективным является метод Ньютона шага 1, выдающий фактически аналитическое решение. На более сложных случаях (особенно для неполиномиальных функций) наилучшим образом себя показала реализация метода Ньютона, основанная на золотом сечении, хотя и стоит заметить, что для подобного результата требовалось большее количество значений функции. В промежуточных вариантах относительная эффективность методов может зависеть от конкретной функции. Кроме того, тесты, не включают в себя рассмотрение реальной вычислительной сложности алгоритмов из-за чего в некоторой степени не были раскрыты преимущества Newton-CG, ведь метод сопряжённых коэффициентов позволяет куда быстрее находить обратную к Гессиану матрицу, пусть и с понижением точности (так же этот плюс по сравнению с нашей реализацией мог раскрыться на существенно многомерных функцией - данные тесты не были включены как раз из-за отсутствия способов объективно замерить время выполнения методов, исключив внешние факторы)
Перейдём к сравнению методов Ньютона с другими методами оптимизации.. Параллельно в этом блоке будет рассматриваться изменение поведения методов в случае вычисления значений приближённо. Для сравнения будут использованы: scipy реализация методов Нелдера-Мида, покоординатный спуск и градиентный спуск на основе золотого сечения, представителем квазиньютоновских методов послужит scipy реализация BFGS Тестирование будет производиться на паре функций - первая рассчитывает приближённое значение производной(вычисленное посредством численных методов), вторая - точное
methods = [
ScipyMethod("Nelder-Mead"),
CoordinateDescent(aprox_dec=maxPrec),
GradientDescent(learning_rate=100, aprox_dec=1e-11),
NewtonBase(1),
Newton(aprox_dec=1e-11, learning_rate=1),
ScipyMethod("Newton-CG"),
ScipyMethod("BFGS")
]
Как и в прошлый раз, начнём рассмотрение с функций тривиального вида
result = Runner.run(methods, [base, base_absolute], point, modules1, precision=defPrec, **TABLE, **NO_VISUALIZE)
LambdaOracul +--------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Method name | StepCount | CallCount | GradientCallCount | HessianCallCount | AbsolutePrecisionCount(1e-07) | AbsolutePrecision | +==================================================+=============+=============+=====================+====================+=================================+=====================+ | ScipyMethod(nelder-mead) | 39 | 145 | 0 | 0 | 38 | 3.89735e-08 | +--------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | CoordinateDescent(300) | 31 | 192 | 0 | 0 | 29 | 2.64895e-08 | +--------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | GradientDescent(100,GoldenRatioMethod,eps=1e-11) | 2 | 130 | 2 | 0 | 1 | 1.51882e-12 | +--------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | NewtonBase(1) | 2 | 0 | 2 | 2 | 1 | 0 | +--------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Newton(1,GoldenRatioMethod,eps=1e-11) | 2 | 110 | 2 | 2 | 1 | 0 | +--------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | ScipyMethod(newton-cg) | 2 | 2 | 2 | 2 | 1 | 8.91646e-13 | +--------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | ScipyMethod(bfgs) | 2 | 5 | 5 | 0 | 2 | 2.67417e-12 | +--------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ SymbolOracul +--------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Method name | StepCount | CallCount | GradientCallCount | HessianCallCount | AbsolutePrecisionCount(1e-07) | AbsolutePrecision | +==================================================+=============+=============+=====================+====================+=================================+=====================+ | ScipyMethod(nelder-mead) | 39 | 145 | 0 | 0 | 38 | 3.89735e-08 | +--------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | CoordinateDescent(300) | 31 | 192 | 0 | 0 | 29 | 2.64895e-08 | +--------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | GradientDescent(100,GoldenRatioMethod,eps=1e-11) | 2 | 130 | 2 | 0 | 1 | 1.53655e-12 | +--------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | NewtonBase(1) | 2 | 0 | 2 | 2 | 1 | 0 | +--------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Newton(1,GoldenRatioMethod,eps=1e-11) | 2 | 110 | 2 | 2 | 1 | 0 | +--------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | ScipyMethod(newton-cg) | 2 | 2 | 2 | 2 | 1 | 0 | +--------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | ScipyMethod(bfgs) | 2 | 5 | 5 | 0 | 2 | 5.68712e-15 | +--------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+
Ожидаемо, все методы сошлись. Методам нулевого порядка потребовалось несоклько больше итераций, нежели другим, однако как заявлять об их вычислительной эффективности, так и неэффективности в данном случае было бы неправильно. Хотя им и потребовалось значительное число вычислений функции, в значениях Гессиана и Градиента они не нуждались. В случаях, когда производную вычислить затруднительно, подобная особенность весьма полезна. Так же стоит обратить внимание на то, что при вычислении производной разностным методом Newton-CG не сошёлся до абсолютной точности. Подобный эффект ожидаем, но продемонстрировать его явно всё же стоило (подобные же эффекты могли появиться и у других методов Ньютона) Кроме того, при вычислении производной приближённо точность BFGS упала на несколько порядков. Так же данный метод требует несколько большего количества вычислений функции и производной. Тем не менее, данная функция слишком проста, чтобы быть показательной и за ростом числа вызовов мы пронаблюдает в последующих функциях. Таким образом, уже наблдаются общие тенденции
point = np.array([102.1, 320.5])
result = Runner.run(methods, [rosenbrok, rosenbrok_absolute], point, modules2, precision=defPrec, **TABLE, **NO_VISUALIZE)
LambdaOracul +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Method name | StepCount | CallCount | GradientCallCount | HessianCallCount | AbsolutePrecisionCount(1e-07) | AbsolutePrecision | +=================================================+=============+=============+=====================+====================+=================================+=====================+ | ScipyMethod(nelder-mead) | 108 | 400 | 0 | 0 | Undefined | 61.3151 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | CoordinateDescent(300) | 499 | 2039 | 0 | 0 | Undefined | 10423.9 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | GradientDescent(10,GoldenRatioMethod,eps=1e-11) | 499 | 29940 | 499 | 0 | Undefined | 314.216 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | NewtonBase(1) | 7 | 0 | 7 | 7 | 6 | 2.49247e-14 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Newton(1,GoldenRatioMethod,eps=1e-11) | 217 | 11935 | 217 | 217 | 216 | 5.26391e-15 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | ScipyMethod(newton-cg) | 75 | 276 | 276 | 149 | 1 | 7.07328e-11 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | ScipyMethod(bfgs) | 195 | 491 | 491 | 0 | 194 | 2.63027e-14 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ SymbolOracul +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Method name | StepCount | CallCount | GradientCallCount | HessianCallCount | AbsolutePrecisionCount(1e-07) | AbsolutePrecision | +=================================================+=============+=============+=====================+====================+=================================+=====================+ | ScipyMethod(nelder-mead) | 108 | 400 | 0 | 0 | Undefined | 61.3151 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | CoordinateDescent(300) | 499 | 2039 | 0 | 0 | Undefined | 10423.9 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | GradientDescent(10,GoldenRatioMethod,eps=1e-11) | 499 | 29940 | 499 | 0 | Undefined | 314.082 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | NewtonBase(1) | 5 | 0 | 5 | 5 | 4 | 1.6782e-14 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Newton(1,GoldenRatioMethod,eps=1e-11) | 217 | 11935 | 217 | 217 | 216 | 1.32075e-14 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | ScipyMethod(newton-cg) | 75 | 276 | 276 | 149 | 1 | 7.07026e-11 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | ScipyMethod(bfgs) | 195 | 495 | 495 | 0 | 194 | 2.77462e-12 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+
Как и следовало ожидать, на функции более сложного вида методы нулевого порядка и Градиентный спуск показали себя не лучшим образом, не сумев сойтись к истинному минимуму. Таким образом, помимо меньшего числа итераций до схождения, плюсом ньтоновских и квазиньютоновских методов стоит признать их более широкуб применимость. Однако стоит обратить внимание на понижение точности BFGS при приближённом вычислении значений. Наиболее ярко влияение приближённых вычислений на результат функции может раскрыть пример ниже
bad_defined_absolute = SymbolOracul(1000 * (x - 100) ** 10 + 100 * (y + 20) ** 10, ['x', 'y'])
bad_defined = LambdaOracul(lambda x, y: 1000 * (np.float64(x) - 100) ** 10 + 100 * (np.float64(y) + 20) ** 10)
min_point = np.array([100, -20])
metrics = [StepCount(), CallCount(), GradientCallCount(), HessianCallCount(),
AbsolutePrecisionCount(0.001, min_point),
MinAbsolutePrecision(min_point)]
conditions = [StepCountCondition(1000), PrecisionCondition(0.00001)]
modules = animations + metrics + conditions
methods2 = [
CoordinateDescent(),
ScipyMethod("Nelder-Mead"),
GradientDescent(learning_rate=10, aprox_dec=0.00001),
NewtonBase(1),
Newton(aprox_dec=1e-14),
ScipyMethod("Newton-CG"),
ScipyMethod("BFGS")
]
oraculs = [bad_defined, bad_defined_absolute]
point = np.array([-100, -200])
result = Runner.run(methods2, oraculs, point, modules, precision=minPrec, **TABLE, **NO_VISUALIZE)
LambdaOracul +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+------------------------+ | Method name | StepCount | CallCount | GradientCallCount | HessianCallCount | AbsolutePrecisionCount(0.001) | MinAbsolutePrecision | +=================================================+=============+=============+=====================+====================+=================================+========================+ | CoordinateDescent(300) | 999 | 113 | 0 | 0 | Undefined | 1 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+------------------------+ | ScipyMethod(nelder-mead) | 35 | 130 | 0 | 0 | 25 | 1.50394e-06 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+------------------------+ | GradientDescent(10,GoldenRatioMethod,eps=1e-05) | 30 | 930 | 30 | 0 | 29 | 1.41434e-06 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+------------------------+ | NewtonBase(1) | 345 | 0 | 345 | 345 | Undefined | 0.0170733 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+------------------------+ | Newton(3,GoldenRatioMethod,eps=1e-14) | 254 | 18288 | 254 | 254 | Undefined | 0.0146697 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+------------------------+ | ScipyMethod(newton-cg) | 35 | 68 | 68 | 68 | Undefined | 0.120232 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+------------------------+ | ScipyMethod(bfgs) | 85 | 183 | 183 | 0 | Undefined | 0.148648 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+------------------------+ SymbolOracul +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+------------------------+ | Method name | StepCount | CallCount | GradientCallCount | HessianCallCount | AbsolutePrecisionCount(0.001) | MinAbsolutePrecision | +=================================================+=============+=============+=====================+====================+=================================+========================+ | CoordinateDescent(300) | 999 | 113 | 0 | 0 | Undefined | 1 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+------------------------+ | ScipyMethod(nelder-mead) | 35 | 130 | 0 | 0 | 25 | 1.50394e-06 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+------------------------+ | GradientDescent(10,GoldenRatioMethod,eps=1e-05) | 30 | 930 | 30 | 0 | 29 | 1.4231e-06 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+------------------------+ | NewtonBase(1) | 128 | 0 | 128 | 128 | 107 | 7.62687e-05 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+------------------------+ | Newton(3,GoldenRatioMethod,eps=1e-14) | 41 | 2952 | 41 | 41 | 31 | 1.62229e-05 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+------------------------+ | ScipyMethod(newton-cg) | 35 | 68 | 68 | 68 | Undefined | 0.120224 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+------------------------+ | ScipyMethod(bfgs) | 95 | 196 | 196 | 0 | Undefined | 0.0984491 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+------------------------+
Методы, не использующие производные не изменили свои показатели, ведь для них . BFGS и NewtonCG не смогли вычислить реальный минимум в случае плохо определённой функции вовсе (стоит заметить, что они они находятся в окрестности минимума, функция весьма плоская там, потому как для Newton-CG, использующему в своей реализации приближение обратного Гессиана, так и для BFGS, исходно ищущего лишь приближение Гессиана, нахожденение минимума подобной функции может быть весьма проблематичным). Интереснее ситуация обстоит с NewtonBase(1) и Newton(3,GoldenRatioMethod,eps=1e-14) - найти минимум с искомой точностью они смогли лишь при точном вычислении значений. Объясняется это тем же образом, что и поведение BFGS/Newton-CG - для некоторых функций крайне важна точность вычислений градиента и Гессиана, потому в случае невозможности предоставления подобных гарантий может иметь смысл воспользоваться методами, не требующих производной/Гессиана
Перейдём, наконец, к не полиномиальным функциям. Конечно, функция МакКормика в тестировании учавствовать не будет - как мы видели ранее, найти её минимум оказалось непосильной задачей даже для наиболее мощных(на текущий момент) наших методов, потому тестировать на ней заведоомо более слабые методы было бы не слишком целесообразно (если быть более точным, данное тестирование было проведено, однако ввиду малой информативности в отчёт включено не было). Чтобы не смотреть вновь на тривиальные случаи, сделаем тестируемую ранее не полиномиальную базовую функцию несколько более экстремальной, сохранив общий её вид
import warnings
warnings.filterwarnings("ignore")
base_non_poly_absolute = SymbolOracul(np.e**((x+0.54719)**2) + np.e**((y+1.54719)**2), ['x', 'y'])
base_non_poly = LambdaOracul(lambda x, y:np.e**(np.float64((x+0.54719)**2)) + np.e**(np.float64((y+1.54719)**2)))
point = np.array([3, 2.5])
abs_metrics6 = metrics_base + [AbsolutePrecisionCount(minPrec, [-0.54719, -1.54719]), AbsolutePrecision([-0.54719, -1.54719])]
modules6 = animations + abs_metrics6 + conditions
methods[2] = GradientDescent(learning_rate=10, aprox_dec=1e-11)
result = Runner.run(methods, [base_non_poly, base_non_poly_absolute], point, modules6, precision=defPrec, **TABLE, **NO_VISUALIZE)
LambdaOracul +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Method name | StepCount | CallCount | GradientCallCount | HessianCallCount | AbsolutePrecisionCount(1e-05) | AbsolutePrecision | +=================================================+=============+=============+=====================+====================+=================================+=====================+ | ScipyMethod(nelder-mead) | 32 | 123 | 0 | 0 | 22 | 2.04878e-08 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | CoordinateDescent(300) | 19 | 144 | 0 | 0 | 18 | 5.98683e-06 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | GradientDescent(10,GoldenRatioMethod,eps=1e-11) | 7 | 420 | 7 | 0 | 6 | 2.02588e-08 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | NewtonBase(1) | 23 | 0 | 23 | 23 | 22 | 4.6004e-15 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Newton(1,GoldenRatioMethod,eps=1e-11) | 23 | 1265 | 23 | 23 | 22 | 6.28822e-15 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | ScipyMethod(newton-cg) | 12 | 23 | 23 | 23 | 1 | 2.34833e-10 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | ScipyMethod(bfgs) | 25 | 53 | 53 | 0 | 24 | 2.22615e-08 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ SymbolOracul +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Method name | StepCount | CallCount | GradientCallCount | HessianCallCount | AbsolutePrecisionCount(1e-05) | AbsolutePrecision | +=================================================+=============+=============+=====================+====================+=================================+=====================+ | ScipyMethod(nelder-mead) | 32 | 123 | 0 | 0 | 22 | 2.04878e-08 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | CoordinateDescent(300) | 19 | 144 | 0 | 0 | 18 | 5.98683e-06 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | GradientDescent(10,GoldenRatioMethod,eps=1e-11) | 7 | 420 | 7 | 0 | 6 | 2.0259e-08 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | NewtonBase(1) | 23 | 0 | 23 | 23 | 22 | 2.22045e-16 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Newton(1,GoldenRatioMethod,eps=1e-11) | 23 | 1265 | 23 | 23 | 22 | 2.22045e-16 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | ScipyMethod(newton-cg) | 12 | 23 | 23 | 23 | 1 | 2.34807e-10 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | ScipyMethod(bfgs) | 25 | 53 | 53 | 0 | 24 | 2.22615e-08 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+
Говоря про эффективность методов, в данном случае, вероятно, наиболее важным является тот факт, что всем удалось сойтись. Что важнее, здесь мы можем увидеть другую крайность: функция настолько экстремальна, что считались ли значения честно или приближённо не имеет значения вовсе. Таким образом, хотя и существуют ситуации, при которых приближение вовсе не окажет влияния на результат или будет незначительным, забывать о возможности ухудшения параметров метода при использовани приблиённых вычислений не стоит. На этом, полагаю, с анализом потенциального влияния приближённых вычислений стоит закончить, ведь основные выводы сделаны, а, как показала уже функция Розенброка, рассчёт точных значений может занимать существенное время. Рассмотрим вторуую не полиномиальную функцию в рамках сравнения эффективности методов
test = LambdaOracul(lambda x, y: np.sqrt(x**(2)+np.sin(x+y)+x**(2)+y**(2)-2 * x *y+((3)/(4))))
test_min = [-0.4442042, -0.6663063]
point = np.array([2200.1, 200.5])
methods3 = [
CoordinateDescent(),
ScipyMethod("Nelder-Mead"),
GradientDescent(learning_rate=10, aprox_dec=0.00001),
NewtonBase(1),
Newton(),
ScipyMethod("Newton-CG"),
ScipyMethod("BFGS")
]
abs_metrics_test = metrics_base + [AbsolutePrecisionCount(minPrec, test_min), AbsolutePrecision(test_min)]
modules_test = animations + abs_metrics_test + [StepCountCondition(500), PrecisionCondition(minPrec)]
result = Runner.run(methods3, [test], point, modules_test, precision=defPrec, **TABLE, **FULL_VISUALIZE)
LambdaOracul +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Method name | StepCount | CallCount | GradientCallCount | HessianCallCount | AbsolutePrecisionCount(1e-05) | AbsolutePrecision | +=================================================+=============+=============+=====================+====================+=================================+=====================+ | CoordinateDescent(300) | 499 | 175 | 0 | 0 | Undefined | 5.88428e-05 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | ScipyMethod(nelder-mead) | 49 | 182 | 0 | 0 | 39 | 2.96122e-08 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | GradientDescent(10,GoldenRatioMethod,eps=1e-05) | 271 | 8401 | 271 | 0 | 270 | 1.96331e-06 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | NewtonBase(1) | 499 | 0 | 499 | 499 | Undefined | 4.27886 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Newton(3,GoldenRatioMethod,eps=0.0001) | 67 | 1608 | 67 | 67 | 66 | 5.11053e-09 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | ScipyMethod(newton-cg) | 10 | 26 | 26 | 19 | 1 | 3.68818e-09 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | ScipyMethod(bfgs) | 16 | 49 | 49 | 0 | 16 | 1.66056e-09 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+
Функция имеет довольно любопытный вид в окрестности своего нуля за счёт добавления синуса. Как можно увидеть, на функции полохо приближаемой квадратичной, метод Ньютона с фиксированным шагом вновь разошёлся. Остальные методы сошлись, хотя стоит отметить аномалию в числе итераций метода Ньютона с переменным шагом. Вероятно, она обусловлена как раз тем, что используется весьма примитивный метод золотого сечения для поиска минимума на луче. Стоит также отметить, что это уже не первый тест, где эффективность выбора направления методом Ньютона была ограничена низкой эффективностью метода выбора точки.
Подводя итоги сравнения и отдельного анализа метдов Ньютона, приходится констатировать, что хотя они и требуют меньшего числа шагов для схождения, их вычислительная эффективность сомнительна. Так, метод BFGS, не требуя Гессиана, выдавал показатели, схожие с методом Newton-CG, кроме того, в случае плохо апрокимируемых функций или работы в области функции, апроксимация которой квадратичной недостаточно точна, методы Ньютона без доверительного региона могут демонстрировать весьма странное поведение. Пример, представленный в дополнительном задании 2 также представляет не с лучшей стороны классичекий метод Ньютона. Многие из этих проблем решаются более сложным алгоритмом выбора точки, тем не менее, подобное усложнение влечёт и увеличение вычислительной сложности, чего хотелось бы избежать (особенно, если речь идёт о вычислениипроизводных). Кроме того, иногда ошибки алгоритма выбора шага могут даже ухудшить работу метода. Один из примеров подобного был представлен ранее, когда метод Newton-CG сойтись не сумел, в то время как это вышло у куда менее интеллектуального NewtonBase. Далее мы подробнее рассмотрим изменения в поведении алгоритма Ньютона, если выполнять одномерный поиск в соответствии с усиленным равилом Вольфе.
Дополнительное задание 1¶
Класс NewtonWolfe, далее именуемый Wolfe реализует метод Ньютона с одномерным поиском по сильному правилу Вольфе. Ранее наши алгоритмы использовали весьма примитивные методы одномерного поиска, что существенно ограничивало их применимость. Однако прежде чем перейти к тестам, приведём описание метода и некоторые особенности реализации
Описание работы правила Вольфе: Как всегда задача минимизировать функцию. Для этого у нас есть гладкая функция f, и мы сведем задачу min(f(x)) -> min(f(xk + a * pk)), pk - направление поиска, а a - длина шага
Есть два вида условий Вольфе: слабое и строгое. Строгое имеет на одино неравенство больше (3 пункт). Перечеслим эти 3 неравенста:
- f(xk + ak * pk) <= f(xk) + c1 * ak * transpose(pk) * grad(f(xk))
- -transpose(pk) * grad(f(xk + ak * pk)) <= -c2 * transpose(pk) * grad(f(xk))
- | transpose(pk) * grad(f(xk + ak * pk)) | <= c2 * | transpose(pk) * grad(f(xk)) |
0 < c1 < c2 < 1. Обычно c1 выбирается довольно маленьнким, а c2 сильно больше. Самые популярные значения c1 = 10e-4, c2 = 0.9. Рассмотрим условия поподробнее: pk выбирается в зависимости от алгоритма поиска. Например, градиентный спуск использует pk равный отрицательному градиенту. В нашем же случае pk выбирается методом Ньютона. В качестве гиперпараметров подаются c1, c2, max_iter, отвечающий за приближение и learning_rate, отвечающий за максимальную длину шага
- Гарантирует, что длина шага уменьшается в достаточной степени
- Гарантирует достаточное уменьшение уклона
- 1 + 3 условия образуют строгое условие Вольфе и они обечспечивают то, что результат будет близок к критической точке
Для внесения разнообразия, в слуае запуска на тех функциях, что уже встречались в тестах, начальные оточки будут изменены. Кроме того, часть функций для большей показательности будет заменена на более примечательные
newton_methods = [
NewtonBase(1),
Newton(learning_rate=4, aprox_dec=defPrec),
ScipyMethod("Newton-CG")
]
all_newton = newton_methods + [NewtonWolfe(c1=0.1, c2=0.9, aprox_dec=1e-7, max_iters=20)]
point = np.array([-100, -2000])
result = Runner.run(all_newton, [rosenbrok], point, modules2, precision=defPrec, **TABLE, **FULL_VISUALIZE)
LambdaOracul +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Method name | StepCount | CallCount | GradientCallCount | HessianCallCount | AbsolutePrecisionCount(1e-07) | AbsolutePrecision | +=======================================+=============+=============+=====================+====================+=================================+=====================+ | NewtonBase(1) | 6 | 0 | 6 | 6 | 5 | 7.28915e-14 | +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Newton(4,GoldenRatioMethod,eps=1e-07) | 124 | 4836 | 124 | 124 | 123 | 2.31022e-15 | +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | ScipyMethod(newton-cg) | 20 | 54 | 54 | 39 | 1 | 3.07729e-10 | +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Wolfe(0.1,0.9,eps=1e-07) | 19 | 28 | 47 | 19 | 18 | 2.45771e-14 | +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+
Начнём серию тестов с довольно простой для методов Ньютона функции - функции Розенброка. Заметим, что по сравнению с другими модификациями метода Ньютона, использующими методы одномерного поиска, основанная на Вольфе сошлась за наименьшее число шагов, да и затратила значительно меньше вычислений. Тем не менее, результат NewtonBase из-за относительной простоты функции обойти не удалось
all_newton = newton_methods + [NewtonWolfe(c1=0.1, c2=0.9, aprox_dec=1e-7, max_iters=20)]
point = np.array([-100, -200])
result = Runner.run(all_newton, [bil], point, modules_bil, precision=minPrec, **TABLE, **FULL_VISUALIZE)
LambdaOracul +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Method name | StepCount | CallCount | GradientCallCount | HessianCallCount | AbsolutePrecisionCount(1e-05) | AbsolutePrecision | +=======================================+=============+=============+=====================+====================+=================================+=====================+ | NewtonBase(1) | 72 | 0 | 72 | 72 | Undefined | 3.04138 | +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Newton(4,GoldenRatioMethod,eps=1e-07) | 499 | 19461 | 499 | 499 | Undefined | 167532 | +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | ScipyMethod(newton-cg) | 6 | 12 | 12 | 11 | Undefined | 90.5308 | +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Wolfe(0.1,0.9,eps=1e-07) | 72 | 73 | 145 | 72 | Undefined | 3.04138 | +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+
print_points(result[0])
('no_show_value', array([8.19025073e-17, 1.00000000e+00]), 'Method name', 'NewtonBase(1)')
('no_show_value', array([-1.67529114e+05, 1.00000592e+00]), 'Method name', 'Newton(4,GoldenRatioMethod,eps=1e-07)')
('no_show_value', array([ 3.56858548e-06, -8.99811040e+01]), 'Method name', 'ScipyMethod(newton-cg)')
('no_show_value', array([8.19025073e-17, 1.00000000e+00]), 'Method name', 'Wolfe(0.1,0.9,eps=1e-07)')
Результаты работы на функции Билла оказались уже несколько более интересными. Метод Вольфе показал себя лучше других методов Ньютона, использующих методы одномерного поиска, сойдясь к точке, находящейся ближе к минимуму. Тем не менее, удивительно, что туда же сошёлся и NewtonBase, притом за то же число итераций
point = np.array([-10, 20.5])
result = Runner.run(all_newton, [hard_non_poly], point, modules4, precision=minPrec, **TABLE, **FULL_VISUALIZE)
LambdaOracul +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Method name | StepCount | CallCount | GradientCallCount | HessianCallCount | AbsolutePrecisionCount(1e-05) | AbsolutePrecision | +=======================================+=============+=============+=====================+====================+=================================+=====================+ | NewtonBase(1) | 186 | 0 | 186 | 186 | Undefined | 2.2518e+15 | +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Newton(4,GoldenRatioMethod,eps=1e-07) | 4 | 156 | 4 | 4 | Undefined | 7.49736 | +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | ScipyMethod(newton-cg) | 3 | 8 | 8 | 5 | Undefined | 7.49736 | +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Wolfe(0.1,0.9,eps=1e-07) | 6 | 10 | 16 | 6 | Undefined | 7.49736 | +---------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+
print_points(result[0])
('no_show_value', array([-1.59226503e+15, 1.59226503e+15]), 'Method name', 'NewtonBase(1)')
('no_show_value', array([5.30143761, 5.30143761]), 'Method name', 'Newton(4,GoldenRatioMethod,eps=1e-07)')
('no_show_value', array([5.3014376, 5.3014376]), 'Method name', 'ScipyMethod(newton-cg)')
('no_show_value', array([5.3014376, 5.3014376]), 'Method name', 'Wolfe(0.1,0.9,eps=1e-07)')
Однако, на примере hard_non_poly мы можем увидеть, что хотя Wolfe и немного теряет эффективность в сравнении с NewtonBase, на сложных функциях он сопособен показать себя значительно лучше. Функция сошлась, пусть и к локальному минимуму
Таким образом, модификация метода Ньютона, основанная на правиле Вольфе действительно выглядит как наиболее эффективной в плане вычислений, так и точной в плане сходимости в сравнении с другими вариациями
Переходим к дополнительному заданию 2¶
Стоит понимать, что метод Ньютона может подразумевать различные реализации, детали которых способны оказывать существенное влияние на работут методов. Так, например, в основном задании была представлена точка на функции Розенброка, при старте их которой методы NewtonBase и Newton, в то время как Newton-CG разошёлся. Потому далее, при аналитических размышлениях, будет рассматриваться метод Ньютона, логика работы которого соответствует представленной на лекции, в испытаниях же будут участвоввать лишь методы NewtonBase и Newton
Разобравшись с терминологией, перейдём к содержательной части. Заметим, что метод Ньютона требует вычисления значения, обратного Гессиану. Данная операция в общем случае может быть невозможна, более того, даже если мы считаем Гессиан ненулевым, его абсолютные значения могут быть весьма малы, что приведёт к излишне большому знаению шага. Фактически, если мы видим, что функция имеет малую кривизну, то приближаем её некорректно, принимая за очень сплюснутую квадратичная. Данную проблему сопсобно в значительной степени решить введение доверительного региона, тем не менее, алгоритмы, использующие данную методику не рассматриваются в рамках данной работы.
Осознав это, нахождение функции, на которой классический алгоритм Ньютона не сумеет сойтись становится не столь проблематчиным: достаточно лишь взять функцию, содержащую область, на которой значение Гессиана крайне мало и точку на этом плато. Конечно, для этого подошла бы и функция Изома из прошлой работы, однако пример не был бы столь показательным, ведь на значительном удалении от минимума у данной функции и значения производной стаовятся бесполезными. Потому рассмотрим следующую функцию:
expr = (sympy.atan(x+2*y) - ((np.pi / 2) - sympy.atan(x - y)))**2
test2_abs = SymbolOracul(expr, ['x', 'y'])
Она имеет две области минимумов, остальные же участки представляют собой слегка наклонённые плато.
test2 = LambdaOracul(lambda x,y: (np.arctanh(x+2*y)-(((np.pi)/(2))-np.arctanh(x-y)))**(2))
point = np.array([-4000, 0])
abs_metrics_test2 = metrics_base
modules_test2 = animations + abs_metrics_test2 + [StepCountCondition(100), PrecisionCondition(minPrec)]
result1 = Runner.run([GradientDescent(aprox_dec=1e-6), NewtonBase(), Newton(aprox_dec=1e-6)], [test2_abs], point, modules_test2, precision=minPrec, **TABLE, **FULL_VISUALIZE)
SymbolOracul +--------------------------------------------------+-------------+-------------+---------------------+--------------------+ | Method name | StepCount | CallCount | GradientCallCount | HessianCallCount | +==================================================+=============+=============+=====================+====================+ | GradientDescent(300,GoldenRatioMethod,eps=1e-06) | 27 | 1161 | 27 | 0 | +--------------------------------------------------+-------------+-------------+---------------------+--------------------+ | NewtonBase(1) | 99 | 0 | 99 | 99 | +--------------------------------------------------+-------------+-------------+---------------------+--------------------+ | Newton(3,GoldenRatioMethod,eps=1e-06) | 99 | 3267 | 99 | 99 | +--------------------------------------------------+-------------+-------------+---------------------+--------------------+
print_points(result1[0])
('no_show_value', array([1821.24653958, 1821.24635699]), 'Method name', 'GradientDescent(300,GoldenRatioMethod,eps=1e-06)')
('no_show_value', array([-1.08422070e+21, -5.04235217e+03]), 'Method name', 'NewtonBase(1)')
('no_show_value', array([-4.00009867e+03, -3.20157090e-19]), 'Method name', 'Newton(3,GoldenRatioMethod,eps=1e-06)')
И вот, на тестах мы видим картину: методам Ньютона сойтись не удалось, в то время как градиентный спуск вполне успешно нашёл область минимума. К сожалению, в данном случае, запуск сначала градиентного спуска не поможет сойтись методу Ньютона, ведь вся функция (кроме областей минимума) - плато, однако в случае, если лишь часть функции представляет собой плато, запуск градиентного спуска для выхода с него и последующего запуска метода Ньютона может быть целесообразным. Так же проблему мог бы решить более интеллектуальный выбор точки на луче. Например, на данной функции, с задачей справляется замена метода золотого сечения на более интеллектуальный метод поиска (при условии выбора за напрпавление антиградиента)
Говоря про второй подпункт второго задания, для его выполнения достаточно взять функцию с огромным числом локальных минимумов, тогда даже малейшее отклонение в выборе окажется существенным. И такая функция у нас есть: функция Леви
def levi_function(x, y):
FACTOR = 5
return (np.sin(3 / FACTOR * np.pi * x)**2
+ (x / FACTOR - 1)**2 * (1 + np.sin(3 * np.pi * y / FACTOR)**2)
+ (y / FACTOR - 1)**2 * (1 + np.sin(2 * np.pi * y / FACTOR)**2))
levi = LambdaOracul(lambda x, y: levi_function(x, y))
abs_metrics5 = metrics_base + [AbsolutePrecisionCount(minPrec, [5, 5]), AbsolutePrecision([5, 5])]
modules5 = animations + abs_metrics5 + [StepCountCondition(500), PrecisionCondition(minPrec)]
point = np.array([102.1, 320.5])
result = Runner.run(methods, [levi], point, modules5, precision=defPrec, **TABLE, **FULL_VISUALIZE)
LambdaOracul +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Method name | StepCount | CallCount | GradientCallCount | HessianCallCount | AbsolutePrecisionCount(1e-05) | AbsolutePrecision | +=================================================+=============+=============+=====================+====================+=================================+=====================+ | ScipyMethod(nelder-mead) | 42 | 169 | 0 | 0 | Undefined | 307.502 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | CoordinateDescent(300) | 27 | 154 | 0 | 0 | 26 | 3.40371e-06 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | GradientDescent(10,GoldenRatioMethod,eps=1e-11) | 499 | 29940 | 499 | 0 | Undefined | 1.65874 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | NewtonBase(1) | 499 | 0 | 499 | 499 | Undefined | 315.437 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | Newton(1,GoldenRatioMethod,eps=1e-11) | 2 | 110 | 2 | 2 | Undefined | 327.344 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | ScipyMethod(newton-cg) | 201 | 410 | 410 | 400 | Undefined | 327.839 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+ | ScipyMethod(bfgs) | 6 | 16 | 16 | 0 | Undefined | 315.002 | +-------------------------------------------------+-------------+-------------+---------------------+--------------------+---------------------------------+---------------------+
Результат соответствует ожиданиям: методы Newton-CG и Newton(1,GoldenRatioMethod,eps=1e-11) оказались в крайне близких точках - и действительно, методы во многом схожи: различия возникают лишь в рассчёте обратного значения Гессиана (однако и производные в данном случае мы берём приближённо) и методе выбора точки на луче (в случае с Newton-CG задействовано правило Армихо, в то время как Newton использует лишь золотое сечение), следующей парой, попавшей однако в разные минимумы являются BFGS и NewtonBase(1). Как причины их расхождения - так и схожести очевидны. Оба метода делают шаг по одному правилу, однако NewtonBase(1) использует в своих вычислениях Гессиан, а BFGS - его приближение (довольно точное, стоит заметить). Нелдер-мид, в силу своей специфичности, ушёл в совершенно иную точку
Неожиданно сошёлся к абсолютному минимуму координатный спуск. Приина очевидна: так уж совпали шаг и начальная точка (а это, смею вас заверить, и правда ни что иное, как совпадение), тем не менее, подобное может послужить напоминанием о том, что иногда и более слабые методы, в силу случая или их специфики, способны выдать лучший результат в конкретном случае.
Координатный спуск шёл по направлению градиента и тоже весьма сильно разошёлся с собратьями, выисляющими или приближающими Гессиан. Тем не менее, в данном случае, из-за слишком высокой специфичности функции данный подход оказался более выгодным (и на самом деле, каждый отдельный перепад весьма хорошо приближается квадратичной функцией, однако попасть в глобальный минимум при таком подходе едва ли удастся